
* Set Directory
cd "...Replication Files"

***************************************************************
********************* ASSEMBLING DATASETS *********************
***************************************************************

*** Assembling dataset for total sector
use "Cleaned_Data.dta", clear

keep totrev OurTotal3 exps cont3 progrev3 invinc3 fundinc3 grprof3 OthInc33 ass_boy ass_eoy liab_boy liab_eoy fundbal perc_cont3 perc_progrev3 perc_invinc3 perc_fundinc3 perc_grprof3 perc_OthInc33 statecode fisyr ntee1 ntmaj12 fips fundbal ass_boy ass_eoy liab_eoy paytax numempestimate HHISample PosRevSample

collapse (sum) totrev OurTotal3 exps cont3 progrev3 invinc3 fundinc3 grprof3 OthInc33 ass_boy ass_eoy liab_boy liab_eoy fundbal perc_cont3 perc_progrev3 perc_invinc3 perc_fundinc3 perc_grprof3 perc_OthInc33 paytax numempestimate, by (fisyr)

save "Aggregated_Total.dta", replace

*** Assembling dataset BY subsector
use "Cleaned_Data.dta", clear

keep totrev OurTotal3 exps cont3 progrev3 invinc3 fundinc3 grprof3 OthInc33 ass_boy ass_eoy liab_boy liab_eoy fundbal perc_cont3 perc_progrev3 perc_invinc3 perc_fundinc3 perc_grprof3 perc_OthInc33 statecode fisyr ntee1 ntmaj12 fips paytax numempestimate HHISample PosRevSample

collapse (sum) totrev OurTotal3 exps cont3 progrev3 invinc3 fundinc3 grprof3 OthInc33 ass_boy ass_eoy liab_boy liab_eoy fundbal perc_cont3 perc_progrev3 perc_invinc3 perc_fundinc3 perc_grprof3 perc_OthInc33 paytax numempestimate, by (fisyr ntmaj12)

save "Aggregated_bySubgroup.dta", replace



***************************************************************
********************* FIGURE 1 ********************************
***************************************************************

use "Cleaned_Data.dta", clear

collapse (mean) perc_cont3 perc_progrev3 perc_invinc3 perc_fundinc3 perc_grprof3 perc_OthInc33 

label var perc_cont3 "Contributions"
label var perc_progrev3 "Program Revenue"
label var perc_invinc3 "Investment Income"
label var perc_fundinc "Special Events"
label var perc_grprof3 "Inventory Gross Profit"
label var perc_OthInc33 "Other"

graph pie perc_cont3 perc_progrev3 perc_fundinc3 perc_OthInc33 perc_invinc3 perc_grprof3 , graphregion(lstyle(none)) plabel(_all percent, size(*1.5) color(white) format(%3.0f)) pie(1, color(gs2)) pie(2, color(gs4)) pie(3, color(gs6)) pie(4, color(gs8)) pie(5, color(gs10)) pie(6, color(gs12))
graph export "Figure 1.pdf", as (pdf) replace


***************************************************************
********************* FIGURE 2 ********************************
***************************************************************

use "Cleaned_Data.dta", clear

keep totrev OurTotal3 exps cont3 progrev3 invinc3 fundinc3 grprof3 OthInc33 ass_boy ass_eoy liab_boy liab_eoy fundbal perc_cont3 perc_progrev3 perc_invinc3 perc_fundinc3 perc_grprof3 perc_OthInc33 statecode fisyr ntee1 ntmaj12 fips fundbal ass_boy ass_eoy liab_eoy paytax numempestimate HHISample PosRevSample

collapse (mean) cont3 progrev3 OurTotal3 perc_cont3 perc_progrev3, by (fisyr)

* Figure 2a
twoway lowess perc_cont3 fisyr, bwidth(0.9) lcolor(gs10) || scatter perc_cont3 fisyr, mcolor(black) xtitle("Year") ytitle("") ///
 graphregion(color(white)) yscale(range(0.1 (0.1) 1)) ylabel(0 "0%" 0.1 "10%" 0.2 "20%" 0.3 "30%" 0.4 "40%" 0.5 "50%" 0.6 "60%" 0.7 "70%" 0.8 "80%" 0.9 "90%" 1.0 "100%") legend(off) ytitle("Percent Revenue from Contributions")
 graph export "Figure 2a.pdf", as (pdf) replace

* Figure 2b 
twoway lowess perc_progrev3 fisyr, bwidth(0.9) lcolor(gs10) || scatter perc_progrev3 fisyr, mcolor(black) xtitle("Year") ytitle("") ///
 graphregion(color(white)) yscale(range(0 (0.1) 1)) ylabel(0(0.1) 1) legend(off) ytitle("Percent Revenue from Program Revenue") ylabel(0 "0%" 0.1 "10%" 0.2 "20%" 0.3 "30%" 0.4 "40%" 0.5 "50%" 0.6 "60%" 0.7 "70%" 0.8 "80%" 0.9 "90%" 1.0 "100%")
 graph export "Figure 2b.pdf", as (pdf) replace


***************************************************************
********************* FIGURE 3 ********************************
***************************************************************

use "Aggregated_Total.dta", clear

* GDP Data from BEA (2012 = 100) https://apps.bea.gov/iTable/iTable.cfm?reqid=19&step=2#reqid=19&step=2&isuri=1&1921=survey
rename fisyr year
merge 1:1 year using "GDP Time Series.dta" 
keep if _merge == 3
drop _merge

gen OurTotal3_adj=OurTotal3 / deflator
gen exps_adj=exps / deflator
gen cont3_adj=cont3 / deflator
gen progrev3_adj=progrev3 / deflator
gen invinc3_adj=invinc3 / deflator
gen fundinc3_adj=fundinc3 / deflator
gen grprof3_adj=grprof3 / deflator

tsset year
tsfilter hp cyclGDP = realgdp
tsfilter hp cyclrev = OurTotal3_adj
tsfilter hp cyclprogrev = progrev3_adj
tsfilter hp cyclcont = cont3_adj
tsfilter hp cyclinvinc = invinc3_adj
tsfilter hp cyclfundinc = fundinc3_adj
tsfilter hp cyclgrprof = grprof3_adj

input str24 revtype 
"Total Revenue"
"Donations"
"Program Revenue" 
"Inventory Gross Profit" 
"Special Events" 
"Investment Income"
end

gen corr = .

pwcorr cyclGDP cyclrev
replace corr = r(rho) if revtype == "Total Revenue"

pwcorr cyclGDP cyclcont
replace corr = r(rho) if revtype == "Donations"

pwcorr cyclGDP cyclprogrev
replace corr = r(rho) if revtype == "Program Revenue"

pwcorr cyclGDP cyclgrprof
replace corr = r(rho) if revtype == "Inventory Gross Profit"

pwcorr cyclGDP cyclfundinc
replace corr = r(rho) if revtype == "Special Events"

pwcorr cyclGDP cyclinvinc
replace corr = r(rho) if revtype == "Investment Income"

separate corr, by(revtype == "Total Revenue")
graph bar corr0 corr1, nofill over(revtype, sort(corr) label(angle(90) labsize(*0.7))) bar(1, color(black)) bar(2, color(gs10)) ///
graphregion(color(white)) ytitle("") legend(off) ytitle("Correlation") b1title("Revenue Type", size(*0.75))
graph export "Figure 3.pdf", as (pdf) replace


***************************************************************
********************* FIGURE 4 ********************************
***************************************************************

use "Aggregated_bySubgroup.dta", clear
sort ntmaj12 fisyr
rename fisyr year

keep OurTotal3 year ntmaj12

reshape wide OurTotal3, i(year) j(ntmaj12) string

merge 1:1 year using  "GDP Time Series.dta" 
keep if _merge == 3
drop _merge

foreach letter in AR BH ED EH EN HE HU IN MU PU RE UN {
replace OurTotal3`letter' = OurTotal3`letter' / deflator
}

tsset year
tsfilter hp cyclGDP = realgdp

input str30 subsector str5 abb
"Arts, Culture, and Humanities" "AR"
"Higher Education" "BH"
"Education" "ED"
"Hospitals" "EH"
"Environment" "EN"
"Health" "HE"
"Human Services" "HU" 
"International" "IN"
"Mutual Benefit" "MU"
"Public and Societal Benefit" "PU" 
"Religion" "RE"
"Unknown" "UN"
"Total" ""
end

gen corr = .

foreach var in AR BH ED EH EN HE HU IN MU PU RE UN {
tsfilter hp cyclrev`var' = OurTotal3`var'
pwcorr cyclGDP cyclrev`var'
replace corr = r(rho) if abb == "`var'"
}

* Correlation for Total Revenue is based on the result in Figure 3
replace corr = -0.0438 if subsector == "Total"

* Create figure
separate corr, by(subsector == "Total")
graph bar corr0 corr1, nofill over(subsector, sort(corr) label(angle(90) labsize(*0.7))) bar(1, color(black)) bar(2, color(gs10)) ///
graphregion(color(white)) ytitle("") legend(off) ytitle("Correlation") b1title("Subsector", size(*0.9))
graph export "Figure 4.pdf", as (pdf) replace


***************************************************************
********************* TABLE 1 ********************************
***************************************************************

use "Cleaned_Data.dta", clear

keep ass_boy OurTotal3 HHISample fisyr statecode ein perc_cont3 perc_progrev3 HHI3 progrev3 cont3

gen year = fisyr
merge m:1 year using  "GDP Time Series.dta"
keep if _merge == 3
drop _merge

merge m:1 statecode year using  "State GDP Data.dta" 
keep if _merge == 3
drop _merge

xtset ein fisyr

gen ass_boy_adj=ass_boy / deflator
gen OurTotal3_adj=OurTotal3 / deflator
gen progrev3_adj=progrev3 / deflator
gen cont3_adj = cont3 / deflator

gen OurTotal3_perc_changeT = 100*2*(OurTotal3_adj - l.OurTotal3_adj)/(abs(OurTotal3_adj) + abs(l.OurTotal3_adj)) 
gen progrev3_perc_changeT = 100*2*(progrev3_adj - l.progrev3_adj)/(abs(progrev3_adj) + abs(l.progrev3_adj)) 
gen cont3_perc_changeT = 100*2*(cont3_adj - l.cont3_adj)/(abs(cont3_adj) + abs(l.cont3_adj)) 

gen stateGDP_cont = stategdp_perc_change * perc_cont3
gen stateGDP_prog = stategdp_perc_change * perc_progrev3
gen stateGDP_HHI = stategdp_perc_change * HHI3

gen GDP_cont = gdp_perc_change * perc_cont3
gen GDP_prog = gdp_perc_change * perc_progrev3
gen GDP_HHI = gdp_perc_change * HHI3

gen GDP_lcont = gdp_perc_change * l.perc_cont3
gen GDP_lprog = gdp_perc_change * l.perc_progrev3
gen GDP_lHHI = gdp_perc_change * l.HHI3

gen GDP_l2cont = gdp_perc_change * l2.perc_cont3
gen GDP_l2prog = gdp_perc_change * l2.perc_progrev3
gen GDP_l2HHI = gdp_perc_change * l2.HHI3

gen GDP_l3cont = gdp_perc_change * l3.perc_cont3
gen GDP_l3prog = gdp_perc_change * l3.perc_progrev3
gen GDP_l3HHI = gdp_perc_change * l3.HHI3

gen ln_ass_boy_adj=ln(ass_boy_adj)

keep if HHISample == 1

* Table 1 - Summary Stats
****************************
xtreg OurTotal3_perc_changeT gdp_perc_change GDP_cont perc_cont3 ln_ass_boy_adj, fe vce(cl ein) 
gen sample=e(sample)

outreg2 using "Table 1" if sample==1, word replace sum(log) keep(OurTotal3_perc_changeT gdp_perc_change stategdp_perc_change perc_cont3 perc_progrev3 HHI3 ln_ass_boy_adj) sortvar(OurTotal3_perc_changeT gdp_perc_change stategdp_perc_change perc_cont3 perc_progrev3 HHI3 ln_ass_boy_adj) auto(2)


* Table 2 - Baseline Results
****************************
xtreg OurTotal3_perc_changeT gdp_perc_change GDP_cont perc_cont3 ln_ass_boy_adj, fe vce(cl ein) 
est store tab2_1

xtreg OurTotal3_perc_changeT gdp_perc_change GDP_prog perc_progrev3 ln_ass_boy_adj, fe vce(cl ein) 
est store tab2_2

xtreg OurTotal3_perc_changeT gdp_perc_change GDP_HHI HHI3 ln_ass_boy_adj, fe vce(cl ein) 
est store tab2_3

outreg2 [tab2_1 tab2_2 tab2_3] using "Table 2", word replace sortvar(gdp_perc_change GDP_cont perc_cont3 GDP_prog perc_progrev3 GDP_HHI HHI3 ln_ass_boy_adj)


* Table 3 - GDP with instruments
****************************
xtivreg OurTotal3_perc_changeT gdp_perc_change ln_ass_boy_adj (GDP_cont perc_cont3 = GDP_lcont l.perc_cont3 GDP_l2cont l2.perc_cont3), fe vce(cl ein) 
est store tab3_1

xtivreg OurTotal3_perc_changeT gdp_perc_change ln_ass_boy_adj (GDP_prog perc_progrev3 = GDP_lprog l.perc_progrev3 GDP_l2prog l2.perc_progrev3), fe vce(cl ein) 
est store tab3_2

xtivreg OurTotal3_perc_changeT gdp_perc_change ln_ass_boy_adj (GDP_HHI HHI3 = GDP_lHHI l.HHI3 GDP_l2HHI l2.HHI3) if HHISample == 1, fe vce(cl ein) 
est store tab3_3

outreg2 [tab3_1 tab3_2 tab3_3] using "Table 3", word replace sortvar(gdp_perc_change GDP_cont perc_cont3 GDP_prog perc_progrev3 GDP_HHI HHI3 ln_ass_boy_adj)


* Table 4 - State GDP
****************************
xtreg OurTotal3_perc_changeT stategdp_perc_change stateGDP_cont perc_cont3 ln_ass_boy_adj, fe vce(cl ein) 
est store tab4_1

xtreg OurTotal3_perc_changeT stategdp_perc_change stateGDP_prog perc_progrev3 ln_ass_boy_adj, fe vce(cl ein) 
est store tab4_2

xtreg OurTotal3_perc_changeT stategdp_perc_change stateGDP_HHI HHI3 ln_ass_boy_adj, fe vce(cl ein) 
est store tab4_3

outreg2 [tab4_1 tab4_2 tab4_3] using "Table 4", word replace sortvar(stategdp_perc_change stateGDP_cont perc_cont3 stateGDP_prog perc_progrev3 stateGDP_HHI HHI3 ln_ass_boy_adj)



***************************************************************
********************* PSR SHARE AND LEVELS ON LEFT HAND SIDE ********************************
***************************************************************

* PSR SHARE
xtreg perc_progrev3 gdp_perc_change, fe vce(cl ein)


* PSR SHARE WITH STATE GDP
xtreg perc_progrev3 stategdp_perc_change, fe vce(cl ein)


* PERC CHANGE IN PSR
xtreg progrev3_perc_changeT gdp_perc_change, fe vce(cl ein)

* PERC CHANG IN CONT
xtreg cont3_perc_changeT gdp_perc_change, fe vce(cl ein)